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Kinematic Simulations of turbulent pair diffusion in planar turbulence with a k~^^^ energy spec- 
trum reproduce the results of the laboratory measurements of Jullien et al Phys. Rev. Lett. 82, 
2872 (1999), in particular the stretched exponential form of the PDF of pair separations and their 
correlation functions. The root mean square separation is found to be strongly dependent on initial 
conditions for very long stretches of times. This dependence is consistent with the topological pic- 
^T*' ture of turbulent pair diffusion where pairs initially close enough travel together for long stretches 

of time and separate violently when they meet straining regions around hyperbolic points. A new 
argument based on the divergence of accelerations is given to support this picture. 
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for increasing values of L/rj, where t is time, e is the kinetic energy rate of dissipation per unit mass and Ga is a 
universal dimensionless constant. 



PACS numbers; 47.27.Eq 47.27.Gs 92.10.Lq 47.27.Qb 

The rate with which pairs of points separate in phase space or in physical space is of central importance to the study 
of dynamical systems. Pairs of points in the phase space of a low-dimensional chaotic dynamical system separate 
exponentially. This is the celebrated butterfly effect: dynamics are extremely sensitive to initial conditions. Pairs 
' of fluid elements in fully choatic flows m also separate exponentially. However, in fully developed homogeneous and 
isotropic turbulence, Richardson's law P| stipulates that fluid element pairs separate on average algebraically and in 
such a way that their separation statistics in a certain range of times are the same irrespective of initial conditions. 
Richardson's law is therefore a remarkable claim of universality. Specifically, it stipulates that in a range of times 

1/2 

where the root mean square separation is larger than the Kolmogorov length-scale rj and smaller than the 

^ . integral length-scale L, is increasingly well approximated by 
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■ Richardson accompanied his empirical law (|lj) with a prediction for the probability density function (PDF) of 

pair separations A. The effective diffusivity approach leading to this prediction was criticised by Batchelor Q who 
developed a different approach leading to (^ but also to a different form of the PDF. Kraichnan derived yet another 
expression for the PDF based on his Lagrangian history direct interaction approximation and so did Shlesinger et al 
on the assumption that turbulent pair diffusion is well described by Levy walks. 

, Setting a{t) = A^ and r = A/ a, the PDFs of A predicted by Richardson Q, Batchelor ||] and Kraichnan Q| 

■ are all of the form 

b : 

. ."^ , P(A,i) - cr-iexp(-ar'3) (2) 

with different values of the dimensionless parameters a and (3. Richardson's prediction for the exponent (3 is f3 = 2/3, 
Batchelor's is /3 = 2 and Kraichnan's is /3 = 4/3. The Lagrangian modelling approach of Shlesinger et al ^ leads to 
a totally different, in fact algebraic, PDF form. More recently Jullien et al |^ reported laboratory measurements of 
P(A, t) which are well fitted by (||) with a « 2.6 and P = 0.5±0.1. These laboratory measurements invalidate the PDF 
predictions of Batchelor, Kraichnan and Shlesinger et al. and might raise a question mark over the PDF prediction 
of Richardson even though they can be considered consistent with it if we account for experimental uncertainties. 
Jullien et al Q also observed that fluid element pairs stay close to each other for a long time until they separate 
quite suddenly, a behaviour which seems qualitatively at odds with the effective diffusivity approach adopted by 
Richardson |^ to derive (||) with P = 2/3. In particular, they measured the Lagrangian autocorrelation function of 
pair separations R{t,T) =< A{t)A{t + t) > for — i < t < and found a Lagrangian pair correlation time Tc « 0.6t 
which is surprisingly long. In this paper we report that Kinematic Simulation (KS) |Q reproduces the experimental 
results of Jullien et al KS is a Lagrangian model of turbulent diffusion which is distinct from Levy walks Q and 
makes no use of Markovianity assumptions so that it cannot be reduced to an effective diffusivity approach such as 
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Richardson's Furthermore, the observation that fluid element pairs travel close together for long stretches of time 
until they separate quite suddenly has in fact already been made using KS . 

KS Lagrangian modelling consists in integrating fluid element trajectories by solving ^^^jp- — u(x(i), t) in synthesised 
velocity fields u(x, t). Statistically homegeneous, isotropic and stationary KS velocity fields are superpositions of ran- 
dom Fourier modes Q . KS velocity fields are gaussian but not delta-correlated in time Q] , and this non-Markovianity 
is an essential ingredient in KS. The Lagrangian measurements of were made in an inverse cascade two-dimensional 
turbulent flow. Our KS velocity field is therefore prescribed to be planar and given by 

m=M 

U = ^ [A„i X k,„ COS (k,„ • X + iOmt) 
m—1 

+Bm X km sin (km • X -I- LUmt)] (3) 

where M — 500 is the number of modes, k^ is a random unit vector (k^ — k^M-m) normal to the plane of the 
flow whilst the vectors and Bm are in that plane. The random choice of directions for the m*'* wavemode is 
independent of the choices of the other wavemodes. Note that the velocity field u is incompressible by construction. 
The amplitudes Am and Bm of the vectors A,„ and are determined by the Kolmogorov energy spectrum E{k) 
via the relations = Bm — E{km)^km where Afcm = [km+i — fcm-i)/2. Finally the unsteadiness frequencies ujm 
are determined by the eddy turnover time of wavemode m, that is uim = 0.5-\/A:j^£'(fcm). 

The Lagrangian measurements in |6j were made when the two-dimensional flow had developed an inverse cascade 
with a well-defined k~^^^ energy spectrum. The energy spectrum we have therefore chosen for this study is E{k) 

^^2/3 k~^/^ in the range ^ < k < ^ and equal to outside this range (u'^ is the total kinetic energy of the 
turbulence). The M = 500 wavenumbers are algebraically distributed between ^ and The eddy turnover time 
at the largest wavenumber — can be considered to correspond to a Kolmogorov time scale r^. 
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FIG. 1. Semi- log plot of ap{r) as a function of r = in the case — 



1691. ^ 



0.1 and ^ 



+ 0.10, X 0.20, * 0.30, 



empty box 0.40, black box 0.5. ^ = 1 and ^f- = 0.10, • 0.20, A 0.30, black triangle 0.4 and v 0.5. The solid line is 

ap[r) ~ e 



The inertial range ratio ^ is 0(10) in the laboratory experiment of 



but here we have also run simulations with 



^ = 10, 100, 1691, 11180, 38748, 250000. For initial pair separations Aq smaller or equal to r] our KS integrations lead 



to ctP(A, t) ~ exp(- 



where r = A/a{t) with 2.6 < a < 3 and 0.46 < P < 0.5 in very good agreement with the 



laboratory results of ||6| and for all the ^ values that we tried (see example in Figure 1). (We record, however, that 
this PDF does seem to depend on the initial separation Aq when Ao > rj.) It has already been noted in |^ that 
KS gives non-gaussian stretched exponential PDFs of pair separations without, however, estimating a and /3. The 
synthetic velocity fields of M lead to the Richardson stretched exponential form with /3 = 2/3. An approach based on 
asymmetric Levy walks [[10[ gives rise to stretched exponential forms of aP{A,t) where (3 can be tuned as a function 
of a persistence parameter. 
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FIG. 2. Lagrangian separation correlation factor R[t,T) / [t) as a function of j in the case 
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2.01, X 1.51, * 1.01, empty box 0.76, black box 0.11, 0.06. ^ = 0.1 and 
black triangle down 0.04. 



= • 1.25, A 0.62, black triangle 0.16, V 0.08, 



Following ||] we also calculate correlation functions of pair separations, i.e. R{t, t) =< A(i) A(t+r) > for — t < r < 
and with Aq equal to rj in one set of runs and O.It^ in another (the choice of Aq in is within this range). The 
laboratory results [D show that R{t, r)/cr^(t) is a function of r/t and exactly the same collapse is found here with KS 
(Figure 2). We calculate a Lagrangian correlation time from R{t,T) in the way done in and we obtain Tc « 0.45t 
from Figure 2. This value 0.45 is indeed the constant asymptotic value that we obtain for all large enough scale ratios 
^7 i-c. ^ > 0(10), and it compares sufficiently well with « 0.6t in the laboratory experiment The agreement 
is therefore good and KS leads to the same conclusion, effectively as the laboratory experiment, that pair separations 
remember about half their history. 
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FIG. 3. Non dimensional diagonal Lagrangian velocity correlation _Dii(t, r)/Dii(i, 0) as a function of r/t. Same case as 
Figure |[ = + 0.03, x 0.06, * 0.08, empty box 0.11, black box 0.25, 0.50, • 0.76, A 1. 
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The last set of statistics measured by [|| are Lagrangian correlations of pair velocity differences, i.e. D 



t) > with — t < T < 0, where V/"{t) denotes the ith component of the Lagrangian relative velocity 
between a pair of fluid elements. We calculate these same statistics using our KS model and find that Dij remains 
close to for i 7^ j, that Dii(t, r)/£)ii(t, 0) and D22{t, t) / D22{t,0) are functions of r/tjsee Figure 3) and that this 
collapse is the same for Dn and D22 again in agreement with the laboratory results of 
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FIG. 4. Pair diffusion as a funtion of time for ^ = 38748 and r,, = 0.0027;;j7 and different initial separations, a) ^^n^'jl^ as 

a function of t^, from top to botom ^ = 1000, 100, 10, 1, 0.1, 0.01 and 0.001. b) as a function of for ^ = 1, 

0.1, 0.01 and 0.001. " " 



Having validated our KS Lagrangian model of pair diffusion in planar turbulence against the laboratory experiment 
of , we now turn our attention to Richardson's law (1) and the claim of universality that it is based on. We do indeed 
observe this law over the entire inertial range of times Trj < t < L/u', but only for initial separations Aq between 
7] and 0.1?7, and this for all the ratios ^ that we tried (see Figure 4a). Of course this ratio should be large enough, 
otherwise Richardson's law is not observed for any Aq, but it is surprising that Richardson's law is so Ao-specific 
even at enormous values of reaching O(IO^). 

When Aq < rj, Richardson's law (1) is observed over the limited large scale range 0.2^ to ^ (Figure 4), and the 
coefficient 0.2 seems to have no dependence on ^ in our simulations as long as L/i] is order 10^ or larger, so if it has 
one it must be weak. In the remainder of the inertial range between r,, and 0.2^, the time dependencies of (A — Ao)^ 

and A^ are different from Richardson's (1) and different for different values of Aq < ry (Figure 4), even at extremely 
high ^. We tried to replace t hy t — to{Ao) , where to(Ao) is a virtual origin significantly smaller than 0-2^ but 
did not recover Richardson's law (1), particularly since the discrepancies we observe are over such wide time ranges. 
When Aq is significantly larger than 77 there is no clear indication of a Richardson law at all (Figure 4). 

We have carefully studied the time dependence of (A — Aq)^ in the range between to for different values of 
Aq < ?7 and have found the following formula to collapse the data in that range (see Figure 4b): 



(A-Ao)2 = Ga— iV(i,Ao) 
where the dimensionless function / is given by (using T = 0.2-^) 



(4) 



f{t, Ao) = exp 



ln[^(Ao^)^ 
21n(T„/r) 



In(VT) - ^Jln'it/T) 



(5) 



Note that f{t, Aq) tends to 1 when t is between T and and L/rj — * 00 (i.e. T/r^ 00). The Richardson constant 
Ga is determined from the value of '}^,jf^°l:i in the range 0.2^ to ^ and we find Ga ~ 0.03 for large enough scale 



(«"/L)ti' 

ratio — (of order 10^ and larger). We should stress that in KS, Ga effectively contains both the original Richardson 

constant as in GAci but also the constant of proportionality relating the kinetic energy dissipation rate to We 
therefore retain the orders of magnitude of Ga obtained by KS but not the actual values. 

Integrations of A^ in Direct Numerical Simulations (DNS) of two-dimensional turbulence in the inverse en ergy 
cascade regime also show a strong dependence on Aq, even when the simulations are very highly resolved 
Nowadays, such DNS cannot reach well-defined —5/3 ranges over more than two decades, i.e. L/j] = O(IOO), and 
this at the very highest resolutions currently available. The deviations from Richardson's law (1) observed when 
L/rj = 0(100) might perhaps be due to edge effects [L/rj too small to reach the asymptotic Richardson's law expected 
to be valid for L/rj ^ 1). But can this also be the case in our KS where L/rj reaches 0(10*^)? Clearly one cannot 
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answer this question with numerical simulations except if in the future KS and/or DNS runs with even higher L/rj 
eventually converge to Richardson's law without Aq - dependencies. 

Nevertheless, the success of our KS to reproduce the laboratory observations of and its failure to retrieve 
Richardson's law without Aq - dependencies even at extremely high L/rj does raise the question of the validity of 
Richardson's universality and of the locality assumption that it is based on even asymptotically for arbitrarily 
high L/rj. In general, A^ is a function of L, ry, Aq and u' in KS, and the Richardson locality assumption adapted for 
KS states that, for large enough L/rj, Jj^^ should only depend on A-^ and E{k) at fc = 2tt / \/ A^ when max{rj, Aq) <C 

\/a^ <^ L. Fung & Vassilicos (1998) found this assumption to be valid in planar KS for different spectral exponents 
p between 1 and 2 {E{k) ^ k^P) but specifically for Aq = rj/2 and unsteadiness parameter A = 0(1) and smaller than 
1 in — Xy/ k^En{k,n) . The direct consequence of this assumption is that A^ ^ V with 7 = which is indeed 
observed in KS for different values of p but only for Aq close to and below 77 . What could invalidate locality and 
Richardson's law for Ag very different from rjl 

The low values of Ga and the very large Lagrangian flatness factors of V/" also observed in KS |Q are consistent 
with the observation that fluid element pairs travel close to each other for long stretches of time and separate in 
sudden bursts Fluid element accelerations a = -^u (where = ^ + u • V) are such that V • a = — ^ 

where s is the strain rate matrix and u> the vorticity vector. Hence, V • a is large and positive most often in straining 
regions around hyperbolic points of the flow where is large and uj'^ close to 0. Close fluid element pairs can separate 
violently where V • a is large and positive, and the separation is effective if the streamline structure of the turbulence 
is persistent enough in time. Hence, such violent separation events will most often occur where close fluid element 
pairs meet hyperbolic points that are persistent enough. 

Based on their KS results which were limited to Aq = ?//2, Fung & Vassilicos (1998) rephrased Richardson's 
locality assumption as follows: "in the inertial range, the dominant contribution to the turbulent diffusivity ^A^ 



comes from straining regions of size V A^; these straining regions are embedded in a fractal-eddy structure of cat's 
eyes within cat's eyes and therefore straining regions exist with a variety of length-scales over the entire inertial range." 
Davila & Vassilicos |lj] have related 7 to the fractal dimension D of this fractal-eddy streamline structure of straining 
regions when Aq is close to and below rj = 4:/D). These results suggest that when Aq is between 77 and O.I77, the 
evolution of fluid element pairs by bursts when they meet straining regions somehow tunes into the straining fractal 
structure of the flow and gives rise to Richardson's law. This requires some persistence of the streamline structure, 
and indeed Richardson's law is lost when the unsteadiness parameter A is made signiflcantly larger than 1 [Q. 

This topological picture of turbulent pair diffusion suggested by results in previous papers and our argument 
concerning V • a could also explain the strong Aq - dependence of A^. As Aq decreases well below rj, the probability 
for fluid element pairs to encounter a hyperbolic point and be separated by it also decreases and can become so small 
for Aq ^ rj that pairs may travel close to each other for very long times. Eventually, at times nearing L/u', the eddy 
turnover time of the turbulence, the two fluid elements will be separated by the unsteadiness of the flow rather than 
by its streamline structure as they will have to become independent at times t 3> L/u'. They therefore largely bypass 
the relatively persistent straining fractal streamline structure of the turbulence and also Richardson's law as a result. 
For initial conditions Aq ^ rj, the argument based on V • a does not apply and the separation of fluid element pairs 
cannot be considered to be dominated by straining events in the vicinity of hyperbolic regions. In the framework of 
the topological turbulent pair diffusion picture, this is consistent with the absence of a Richardson law for Aq » rj. 
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